Nontrivial Velocity Distributions in Inelastic Gases 



O 

o 

(N 

> 
O 



o 



C/3 



I 

O 

o 



> 
O 



I 

o 
o 



P. L. Krapivsky^ and E. Ben-Naim^ 
^ Center for Polymer Studies and Department of Physics, Boston University, Boston, MA 02215 
^ Theoretical Division and Center for Nonlinear Studies, Los Alamos National Laboratory, Los Alamos, NM 87545 

We study freely evolving and forced inelastic gases using the Boltzmann equation. We consider 
uniform collision rates and obtain analytical results valid for arbitrary spatial dimension d and 
arbitrary dissipation coefRcient e. In the freely evolving case, we find that the velocity distribution 
decays algebraically, P{v,t) ~ íj"'' for sufRciently large velocities. We derive the exponent a{d,e), 
which exhibits nontrivial dependence on both d and e, exactly. In the forced case, the velocity 
distribution approaches a steady-state with a Gaussian large velocity tail. 
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The dynamics of granular gases present novel chal- 
lenges, previously not encountered in ñuid dynamics 
0. Specifically, the strong underlying energy dissipa- 
tion leads to clustering instabilities and strong velocity 
correlations A series of recent experimental and 

theoretical studies reveals a rich phenomenology. In par- 
ticular, velocities are characterized by anomalous statis- 
tics, sensitive to the details of the driving conditions, thc 
density, and the degree of dissipation ■ 

Kinetic theory provides a systematic framework for de- 
riving macroscopic properties from the microscopic coUi- 
sion dynamics [p^p^. Yet, analysis of the corresponding 
Boltzmann equation often involves uncontroUed approxi- 
mations or use of nearly Maxwellian distributions. Moti- 
vated by the latter issue, we examine both unforced and 
forced inelastic gases using a simplified Boltzmann equa- 
tion. Specifically, we shall employ Maxwell's coUision 
rate which is proportional to the typical velocity rathcr 
than the relative velocity , and show that in general, 
it is analyticaUy tractable. 

Previous studies reveal that this kinetic theory leads 
to interesting behaviors in the freely evolving case. In 
one dimension, while moments of the velocity distribu- 
tion exhibit muUiscaUng [ pO[ , the velocity distribution 
itself stiU approaches a scaling form with an algebraic 
large velocity tail [^. An algebraic tail was also found 
numerically in two dimensions pl¡. Here, we show ana- 



lytically that in arbitrary spatial dimension the velocity 
distribution admits a scaling solution with an algebraic 
large velocity tail. The corresponding exponent, a root 
of a transcendental equation, depends on the spatial di- 
mension and the restitution coefficient. In the driven 
case, we find that although the velocity distribution is 
non-MaxwelUan, it does exhibit a Gaussian tail. 

Our starting point is a homogeneous system of identi- 
cal inelastic spherical particles in arbitrary spatial dimen- 
sion d. The mass and the cross-section are set to unity 
without loss of generality. When two particles collide, 
the normal component of the relative velocity is reduced 
by the restitution coefficient r = 1 — 2e, while the tan- 
gential component remains the same. Denoting by n the 
unit vector connecting the centers of the colUding parti- 
cles, the post-collision velocities vi_2 are given by a linear 
combination of the precolUsion velocities Ui^^j 



vi,2 = ui,2 T (1 - e) (g • n) n. 



(1) 



with the relative velocity g = Ui — u^ . The energy dissi- 
pated in a coUision equals AE = — e(l — e)(g • n)^. When 
e = colUsions are elastic, while when e = 1/2 colUsions 
are perfectly inelastic with maximal energy dissipation. 

Let P(v, t) be the normalized density of particles with 
velocity v at time t. Assuming molecular chaos or perfect 
mixing, we arrive at the foUowing Boltzmann equation for 
tUe velocity distribution 



dP{v,t) 
dt 



= (g) J dnj dui J du2 P(ui,í)P(u2,í){<5[v-Ui + (l-e)(g-n)n] -á(v-Ui)j. 



(2) 



X 



The coUision rate in the MaxweU approximation repre- 
sent the typical velocity scale, (g) = \¡T, with T = ^(í^^) 
the "granular temperature" or the average kinetic energy 
per degree of freedom (for hard spheres, the colUsion rate 
equals the actual relative velocity g • n). This evolution 
equation naturaUy describes a stochastic process where 
randomly chosen pairs of particles undergo inelastic col- 
Usions according to (|]) with a randomly chosen impact 
direction n. In writing Eq. (|^) we tacitly ignored the re- 
striction g • n > on the integration range, because the 
integrand obeys the refiection symmctry n n. Thc 



integration measure should be normaUzed, J dn = 1. 

In the absence of energy input the system "cools" 
indefinitely according to Haff's law [|2^. Indccd, the 
time dependence of the temperature is found from the 
Boltzmann equation (||) to give = -AT^/^ with 

A = 2e(l — e) J dnni (the first axis was conveniently cho- 
sen to be parallel to g). Since nf + . . . + n'^ = 1 and 
J dn = 1 one has A = 2e(l — e)/d. Thus, the tempera- 
ture decays according to T(í) = Tq (1 + í/íq)"^, with the 
initial temperature Tq and the characteristic time scale 
íq = rf/[e(l — e)VTo] • The temperature quantifies veloc- 



1 



ity ñuctuations. In the following, we focus on the natural 
case of isotropic velocity distributions. We wiU show that 
asymptotically, the temperature represents the only rele- 
vant velocity scale as the velocity distribution approaches 
the scaling form 



P(v,í) 



1 



rpd¡2' 



(3) 



Given the convolution structure of the Boltz- 
mann equation (H), we introduce P(k, í), the Fourier 
transform of the velocity distribution function, 
F(k, t) = J dv e*'' '^ P(v, í). Applying the Fourier trans- 
form to Eq. (0) and integrating over the velocities gives 



-^^P(k,í)-HP(k,í) = j dnP[k-q,í]P[q,í], 



(4) 



with q — (1 — e)(k • n)n. This rate equation reflects 
the momentum transferred between particles during col- 
lisions. We seek an isotropic scaling solution for the 
Fourier transform, thc equivalent of (0), 



F(k,í) = $ {k^T) 



(5) 



In the /c — í- limit, the Fourier transform, 
i T, implies that the first two terms in 



with = |k 
F(k,í)-1 

the Taylor expansion of the corresponding scaling func- 
tion ^{x) are universal, ^{x) = 1 — ^x. 

Let us first consider the simpler ID case. Equation (^) 
reduces to -^-§^F{k,t) + F{k,t) = F[ek,t] F[k - €k,t] 
and the scaling function (^) satisfies [Q 

- Xx<P'{x) + <í>{x) = <P[e^x] <Í>[{1- e)'^x] . (6) 

This equation admits a very simple solution [El[ 



$(a;) = [1 + V^] e-^. 



(7) 



The inverse Fourier transform gives the scaling function 
of the velocity distribution as a squared Lorentzian 



V{w) 



1 



TT {l+W^)"^' 



(8) 



Therefore, thc form of the scaling solution (^ is universal 
as it is independent of the dissipation degree. Another 
important feature is the algebraic tail of the velocity dis- 
tribution, V{w) ~ w~^ as ií; — s- oo. 

We now return to the general d-dimensional case. Sub- 
stituting (B|) into (^) and using the temperature cooling 
equation §T = -AT^/^ one flnds that the scaling func- 
tion ^{x) satisfies $(x) - Xx^'{x) = Jdn^{^x)^{T]x), 
where ^ = 1 — (1 — e'^)nf and r; = (1 — e)'^n\. To 
integrate over the impact direction n we use spher- 
ical coordinates and treat the flrst axis as the po- 
lar axis, ni = cos^. The é'-dependent factor of 
the measure is dn = N""^ {sm 9)'^-'^d9 with the factor 
N = {sin 6)'^-^ d6 = B {^ , ensuring proper nor- 

malization {B{a,b) is the beta function). Using ¡j, = 



cos^ 9 as the integration variable, the above governing 
equation for the scaling function reads 



Xx<S>'{x) + <^{x) = / Vfi <S>{£,x) <^{r]x). 



(9) 



where ^(/i) = 1 — (1 — e^)/x and r]{fi) = (1 — e)^/i. Ad- 
ditionally, the integration measure was re-written using 
Vfi, defined via B (i, V]i = fi^^ (1 — fi)^~ dfj. (it 

remains normalized, J^ Vfi = 1). In thc elastic case, the 
velocity distribution is Maxwellian. Indeed, ^ + rj = 1 
and A = when e = 0, and thence ^{x) = e~^/^. 

Our primarily goal is to determine statistics of ex- 
tremely fast particles, namely the tail of the velocity 
distribution. This can be accomplished by noting that 
the large-i" behavior of the velocity distribution is re- 
flected by the small-fc behavior of its Fourier transform. 
This is seen from the small-a: expansion of the ID solu- 
tion (|^) which contains both regular and singular terms: 
<i>(x) = 1 — ix+i x^^^ + • • •. The dominant singular 
x^/^ term reflects the w;~* tail of V{w). In general, an 
algebraic tail of the velocity distribution (0) , 



V{w) 



when 



oo. 



(10) 



indicates the existence of a singular component in the 
Fourier transform, 



$sing(a;) - x'^"'-'^^/^ when 



a; ^ 0, 



(11) 



and vice versa. To see this, it is convenient to take 
the Fourier transform <i>(a;) cx J^ dw w'^^^V^w) e™^ , 
and write x = —s^ thus recasting it into the Laplace 
transform I{s) oc J^ dw w'''~^V{w) e~^'' . The small-s 
expansion of the integral I{s) contains regular and sin- 
gular components. For example, when a < d, the inte- 
gral I{s) diverges as s — > and integration over large-w 
yields the dominant contribution /sing(s) ^ s'^~'^. When 
d < (7 < d+1, 1{0) is flnite, but the next term is the above 
singular term, so I{s) = I{0) + /sing(s) 
the singular contribution is IsinAs) ^ s° 
ing to the singular term of Eq. ( |ll| ) . This is of course con- 
sistent with the one-dimensional behavior where cr = 4. 

The exponent a can now be obtained by inserting 
^{x) = <i>rcg(a;) + $sing(2;) into Eq. (^ and equating the 
dominant singular terms. Combining <i>rcg(0) = 1, with 
thc anticipated leading singular term of Eq. (nW , we flnd 
that the exponent ct is a root of the foUowing integral 
equation 



In general, 
''^, thereby lead- 



1-A^= / Vf, 



^{a-d)/2 ^^{a-d)/2 



(12) 



This equation has a trivial solution a = d + 2, foUow- 
ing from the identity J Vfi {S. + rf) = 1 — A, where the 
singular and the regular components simply coincide, 
j;(<^-d.)/2 _ Since we seek the leading singular term, 
the solution of Eq. ([12) must therefore satisfy a > d + 2. 



2 



The integral equation can be written explicitly in 
terms of special functions 



l-e(l-e 
2F1 



-d 



(13) 



d — a 1 d 2 



-d+l'- 



r(f)r(i] 



with 2Fi{a,b; c; z) the hypergeometric function |g3[. In- 
terestingly, the exponent cr = a{d, e) is a root of the 
transcendental equation (13) and thence it depends in a 
nontrivial fashion on spatial dimension d and the resti- 
tution coefñcient e. 

We first consider the dependence on the restitution co- 
efñcient by considering the quasi-elastic hmit e ^ 0. As 
discussed above, in the elastic case the velocity distribu- 
tion is MaxweUian and the Fourier transform is simply 
$(a;) = e~^/^ |^^. This imphes a diverging exponent 
cr — !■ CX3 as e 0. Therefore, the right-hand side of 
Eq. (O) vanishes, and the leading behavior is 



d 

(7 ~ — 

e 



as 



0. 



(14) 



One can further expand a{d, e) in the e — !■ hmit to 
find: cr(d,e) = de'^ + ai{d) e~^/^ + a^^d) + ■ ■ ■. We 
merely quote the leading correction in the physically rel- 
evant spatial dimensions ai(2) = — 2(e~^ -I- 1)/^/?? and 
ai(3) = —^JZ-K ¡2. Clearly, the quasi-elastic hmit is sin- 
gular. Dissipation, even if minute, seriously changes the 
nature of the system |^,|lO|,^ . 

Next, we discuss the depcndence on the dimension. 
First, one can verify that cr = 4 whcn d = \ using the 
identity 2Fi{a,b;b; z) = (1 — z)^"". The case d = 1 is 
unique in that the entire scahng function and in partic- 
ular the exponent become independent of the restitution 
coefficient. The case of cZ ^ 00 is similar to the e = 
case in that the inelastic nature of the colfisions becomes 
irrelevant, the velocity distribution is MaxweUian, and 
the exponent diverges, cr — > 00 as c¿ — > 00. In this limit, 
the second integral in Eq. ( [l^ is neghgible as it van- 
ishes exponentiaUy with the dimension. The first inte- 
gral can be evaluated by taking the limits d — > 00 and 
/i — > 0, with z = /id/2 being finite. Then, the integra- 



of the energy dissipation rate A = 2e(l — e)/d on d 
and e. Hence, the completely inelastic case provides a 
lower bound for the exponent, cr(d, e) > cr(d, e = 1/2) 
with a{d, 1/2) = 6.28753, 8.32937, for d = 2, 3, re- 
spectively. The former value should be compared with 
cr w 5, obtained from numcrical simulations |^l| . The 
algebraic tails are characterized by unusually large ex- 
ponents which may be difficuU to measure accurately in 
practice (for typical granular particles e ~ 0.1 yielding 
a « 30). Figure 1 also shows that the quantity a/d 
weakly depends upon the dimension, and the large-cZ hmit 
( ITsl ) provides a good approximation cvcn for smaU dimen- 
sions. 




FIG. 1. The exact exponent a, obtained from Eq. (|l^, ver- 
sus the dissipation parameter e. The exponent was scaled by 
the dimension d. Shown also is the limiting large dimension 
expression (|l^). 

Thus far, we discussed only freely coohng systems 
where the energy decreases indefinitely. In typical exper- 
imental situations, however, the system is suppUed with 
energy to balance the energy dissipation pO|-p^. The- 
oreticaUy, it is natural to consider white noise forcing 
i.e., couphng to a thermal heat bath which leads to a 
nonequiUbrium steady state. Specifically, we assume that 
in addition to changes due to coUisions, velocities may 



tion measure is transformed P/i ^ {ttz) ^l^e ^ dz, and ^jg^ ^^^^^^ external forcing: ^ |heat = with 



A?x = /o°°dz(7rz)~i/2e-[i+(i 
Performing the integration yields 



Eq. (|12| ) becomes 1 
where u = 5 ~ 1- 

d 

1 — Atí = [1 + (1 — e^)w]~^/^, from which we find u and 



1 + |g - g3 _ ^1/2 (1 ^ 

e(l-e2) 



1/2 



(15) 



j = 1, . . . , d. We use standard uncorrelated white noise 
{^i{t)^j{t')) = 2DSijS{t — t') with a zero average {^j) = 
0. The temperature rate equation is modified by the ad- 
ditional source term -^T + \T^/^ = D, and the system 
approaches a steady state, T^ = {D/X)-^/^ . The relax- 



-T 



,-t/r 



as d ^ cxD. In general, a (x d, and therefore, the alge- 
braic decay becomes sharper as the dimension increases. 
The exponent a{d, e) increases monotonically with in- 
creasing dimension, and additionaUy, it increases mono- 
tonicaUy with decreasing e (see Fig. 1). Both features 
are intuitive as they mirror the monotonic dependence 



ation toward this state is exponential, |To, 

Uncorrelated white noise forcing amounts to diffusion 
in velocity space, and Eq. (Q) is modified as foUows 
In the steady state, the Fourier 



1 d 
Vrdt 

transform, T(k, í = 00) = í'(y) with y = Dk'^, obeys 



VT dt 



{l + y)-^{y) = {^{^y)-^{m)), 



(16) 



3 



where integration with respect to the measure V^, is de- 
noted by (/) = /g 

Equation (^6|) is solved recursively by employing the 

cumulant expansion í'(y) = exp X]n>i(^y)"^i 



ing \ + y — exp 



we recast Eq. 



exp 



E 

n=l 



(-y)" 



= ( exp 



-^(-y)"G„ 



n=l 



Writ- 
into 

(17) 



where G„ = F„(l — — 77"). The desired cumulants F„ 
are obtained by solving for (G^) recursively and then us- 
ing Fn = (G„) /(1 — — 'y")- oí^e dimcnsion (/i") = 1 
and onc immediately obtains (G„) = [^. In higher 
dimcnsions, the averages acquire non-trivial dependence 
on n though the quahtative nature of the solution re- 
mains the same since (G„) and hence the cumulants 
(as well as the moments of the distribution) are positive. 
This imphes the foUowing small-fc behavior of thc Fourier 
transform at the steady state: Ini^oo(k) ~ — (fc/A;o)^. 
Consequently, the large velocity tail of the velocity dis- 
tribution is Gaussian 



Poo{v) 



-(v/vor 



when w 3> wq- 



(18) 



This behavior agrees with the prediction of Ref. |2^] 
derived via a small-e expansion, but diffcrs from the 
stretched exponential behavior exp(— w'^/^) found for the 
driven inelastic hard sphere gas . 

We have studied inelastic gases within the framework 
of the Boltzmann equation with a uniform coUision ker- 
nel. In the freely evolving case, we have shown analyt- 
icaUy that the density of high-energy particles is sup- 
pressed algebraicaUy. The algebraic tails are character- 
ized by remarkably large exponents, and may be hard 
to distinguish from (stretched) exponential tails. Our 
results, combined with previous kinetic theory studies 
which find exponential, stretched exponential, and Gaus- 
sian tails, indicate that the extremal characteristics can 
be very sensitive to parameters such as the restitution co- 
efñcient, and the dimcnsion [^,^,^. These findings are 
aU the more intriguing when compared with molecular 
dynamics simulations and scahng studies which suggest 
universal velocity distributions, at least in some cases 
On the other hand, our results in the forced case 
support the near-MaxweUian assumptions typicaUy used 
to obtain macroscopic transport coefficients from kinetic 
theory [0. 
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